Lorenz attractor - simulating a multidimensional dynamic system

Lorenz Attractor is a famous multidimensional dynamic system often used as an example for chaotic system, i.e a system in which small deviation in initial states expand and amplified which bring larger changes toward the final outcome.

Mathematically, the lorenz attractor is expressed as a set of 3 differential equations as follow:

\(\displaystyle \frac{dx}{dt} = \sigma(y - x)\)

\(\displaystyle \frac{dy}{dt} = x(\rho- z) - y\)

\(\displaystyle \frac{dz}{dt} = xy - \beta z\)

First, let’s define the function with the above formulas and load the necessary packages for solving numerical solutions for these equations.

library(deSolve)
library(plotly)
library(magrittr)

lorenz <- function(t, state, parameters) {
  with(as.list(c(state, parameters)), {
    dx <- sigma * (y - x)
    dy <- x * (rho - z) - y
    dz <- x * y - beta * z
    return(list(c(dx, dy, dz)))
  })
}

To solve these DEs in series of numerical solutions, we will define the initial state variables x, y, z as well as the values of sigma, rho and beta constants, then call the ODE function to provide solutions for a given time series.

# set the initial state and parameters for the system
initial <- c(x = 1, y = 1, z = 1)
parameters <- c(sigma = 10, beta = 8/3, rho = 28)
times <- seq(0, 100, by = 0.01)

# calling ODE solver and save the result
result <- as.data.frame(ode(y = initial, times = times,
                            func = lorenz, parms = parameters))

We can use plotly package to visualize the lorenz attractor interactively, which nicely takes the shape of the famously inspired name, the Butterfly Effect.

plot_ly(data = result, x = ~x, y = ~y, z = ~z) %>%
  add_paths(color = ~time) %>% 
  layout(scene = list(camera = list(eye = list(x = -1, y = 1, z = 0))))